library(ggplot2)
library(ggpubr)
library(tidyverse)
library(broom)
library(AICcmodavg)
library(rstatix)
library(apa)
library(ez)
library(afex)
library(dplyr)


df <- read.csv("data_with_media_measures.csv")
df$partisanship3[df$partisanship2 < 4] <- "Democrat"
df$partisanship3[df$partisanship2 == 4] <- "Independent"
df$partisanship3[df$partisanship2 > 4] <- "Republican"

df <- df[!is.na(df$partisanship3),]
df$partisanship3 <- factor(df$partisanship3, levels = c("Republican", "Independent", "Democrat"))

## standardized
df$partisan_news_flow <- (df$partisan_news_flow - 
                                mean(df$partisan_news_flow, na.rm=TRUE))/sd(df$partisan_news_flow, na.rm=TRUE)

## standardized
df$top_media_slant <- (df$top_media_slant - 
                            mean(df$top_media_slant, na.rm=TRUE))/sd(df$top_media_slant, na.rm=TRUE)


df$far_left <- (df$far_left - 
                        mean(df$far_left, na.rm=TRUE))/
  sd(df$far_left, na.rm=TRUE)
df$far_right <- (df$far_right - 
                        mean(df$far_right, na.rm=TRUE))/
  sd(df$far_right, na.rm=TRUE)

options(scipen = 100)
#df1 <- df[df$partisanship3 != "Independent", ]


################## ANOVA

# overall slant
# boxplot
ggplot(df) +
  aes(x = partisanship3, y = partisan_news_flow) +
  geom_boxplot() +
  xlab("Partisan group") + 
  ylab("News use slant (left to right)") +
  theme_bw()

# descriptive
group_by(df, partisanship3) %>%
  summarise(
    mean = mean(partisan_news_flow, na.rm = TRUE),
    sd = sd(partisan_news_flow, na.rm = TRUE)
  )


# anova
aov1 <- aov(partisan_news_flow ~ partisanship3,
               data = df
)
summary(aov1)

TukeyHSD(aov1)
plot(TukeyHSD(aov1))

post_test <- glht(aov1,
                  linfct = mcp(partisanship3 = "Tukey"))

summary(post_test)

# outlet diversity 
# boxplot
ggplot(df) +
  aes(x = partisanship3, y = outlet_diversity) +
  geom_boxplot() +
  xlab("Partisan group") + 
  ylab("News use diversity") +
  theme_bw()

# descriptive
summary(df$outlet_diversity)
group_by(df, partisanship3) %>%
  summarise(
    mean = mean(outlet_diversity, na.rm = TRUE),
    sd = sd(outlet_diversity, na.rm = TRUE)
  )


# anova
aov2 <- aov(outlet_diversity ~ partisanship3,
            data = df
)
summary(aov2)

TukeyHSD(aov2)
plot(TukeyHSD(aov2))


post_test <- glht(aov2,
                  linfct = mcp(partisanship3 = "Tukey"))

summary(post_test)


# right extremity
# boxplot
ggplot(df) +
  aes(x = partisanship3, y = far_right) +
  geom_boxplot() +
  xlab("Partisan group") + 
  ylab("Right extremity") +
  theme_bw()

# descriptive
summary(df$far_right)
group_by(df, partisanship3) %>%
  summarise(
    mean = mean(far_right, na.rm = TRUE),
    sd = sd(far_right, na.rm = TRUE)
  )


# anova
aov3 <- aov(far_right ~ partisanship3,
            data = df
)
summary(aov3)

TukeyHSD(aov3)
plot(TukeyHSD(aov3))
post_test <- glht(aov3,
                  linfct = mcp(partisanship3 = "Tukey"))

summary(post_test)


# left extremity
# boxplot
ggplot(df) +
  aes(x = partisanship3, y = far_left) +
  geom_boxplot() +
  xlab("Partisan group") + 
  ylab("Right extremity") +
  theme_bw()

# descriptive
summary(df$far_left)
group_by(df, partisanship3) %>%
  summarise(
    mean = mean(far_left, na.rm = TRUE),
    sd = sd(far_left, na.rm = TRUE)
  )


# anova
aov4 <- aov(far_left ~ partisanship3,
            data = df
)
summary(aov4)

TukeyHSD(aov4)
plot(TukeyHSD(aov4))
post_test <- glht(aov4,
                  linfct = mcp(partisanship3 = "Tukey"))

summary(post_test)

##### appendix tests: top slant
# boxplot
ggplot(df) +
  aes(x = partisanship3, y = top_media_slant) +
  geom_boxplot() +
  xlab("Partisan group") + 
  ylab("Right extremity") +
  theme_bw()

# descriptive
summary(df$top_media_slant)
group_by(df, partisanship3) %>%
  summarise(
    mean = mean(top_media_slant, na.rm = TRUE),
    sd = sd(top_media_slant, na.rm = TRUE)
  )


# anova
aov1.2 <- aov(top_media_slant ~ partisanship3,
            data = df
)
TukeyHSD(aov1.2)

##### appendix tests: slant diversity
# boxplot
ggplot(df) +
  aes(x = partisanship3, y = type_diversity) +
  geom_boxplot() +
  xlab("Partisan group") + 
  ylab("Right extremity") +
  theme_bw()

# descriptive
summary(df$type_diversity)
group_by(df, partisanship3) %>%
  summarise(
    mean = mean(type_diversity, na.rm = TRUE),
    sd = sd(type_diversity, na.rm = TRUE)
  )


# anova
aov2.2 <- aov(type_diversity ~ partisanship3,
              data = df
)
TukeyHSD(aov2.2)



# ttest
# overall slant

res <- t.test(partisan_news_flow ~ partisanship3, data = df1, var.equal = TRUE)
res

t_apa(
  res,
  es = "cohens_d",
  es_ci = FALSE,
  format = c("text"),
  info = FALSE,
  print = TRUE
)

mean(df1$partisan_news_flow[df1$partisanship3 == "Republican"])
mean(df1$partisan_news_flow[df1$partisanship3 == "Democrat"])
sd(df1$partisan_news_flow[df1$partisanship3 == "Republican"])
sd(df1$partisan_news_flow[df1$partisanship3 == "Democrat"])
df1 %>% cohens_d(partisan_news_flow ~ partisanship3, var.equal = TRUE)


#abs
res <- t.test(abs(partisan_news_flow) ~ partisanship3, data = df1, var.equal = TRUE)
res
mean(abs(df1$partisan_news_flow[df1$partisanship3 == "Republican"]))
sd(abs(df1$partisan_news_flow[df1$partisanship3 == "Republican"]))
mean(abs(df1$partisan_news_flow[df1$partisanship3 == "Democrat"]))
sd(abs(df1$partisan_news_flow[df1$partisanship3 == "Democrat"]))
df1 %>% cohens_d(partisan_news_flow ~ partisanship3, var.equal = TRUE)

sub1 <- df1[which((df1$partisanship3 == "Republican" & df1$top_media_slant>=0)|
                    (df1$partisanship3 == "Democrat" & df1$top_media_slant<0)),]
nrow(sub1)
res <- t.test(abs(partisan_news_flow) ~ partisanship3, data = sub1, var.equal = TRUE)
res
mean(abs(sub1$partisan_news_flow[sub1$partisanship3 == "Republican"]))
sd(abs(sub1$partisan_news_flow[sub1$partisanship3 == "Republican"]))
mean(abs(sub1$partisan_news_flow[sub1$partisanship3 == "Democrat"]))
sd(abs(sub1$partisan_news_flow[sub1$partisanship3 == "Democrat"]))

t_apa(
  res,
  es = "cohens_d",
  es_ci = FALSE,
  format = c("text"),
  info = FALSE,
  print = TRUE
)


# top slant

## abs
res <- t.test(abs(top_media_slant) ~ partisanship3, data = sub1, var.equal = TRUE)
res
mean(abs(sub1$top_media_slant[sub1$partisanship3 == "Republican"]))
sd(abs(sub1$top_media_slant[sub1$partisanship3 == "Republican"]))
mean(abs(sub1$top_media_slant[sub1$partisanship3 == "Democrat"]))
sd(abs(sub1$top_media_slant[sub1$partisanship3 == "Democrat"]))

t_apa(
  res,
  es = "cohens_d",
  es_ci = FALSE,
  format = c("text"),
  info = FALSE,
  print = TRUE
)
## all
res <- t.test(top_media_slant ~ partisanship3, data = df1, var.equal = TRUE)
res
t_apa(
  res,
  es = "cohens_d",
  es_ci = FALSE,
  format = c("text"),
  info = FALSE,
  print = TRUE
)

sd(df1$top_media_slant[df1$partisanship3 == "Republican"])
sd(df1$top_media_slant[df1$partisanship3 == "Democrat"])
df1 %>% cohens_d(top_media_slant ~ partisanship3, var.equal = TRUE)


# outlet diversity
res <- t.test(outlet_diversity ~ partisanship3, data = df1, var.equal = TRUE)
res
sd(df1$outlet_diversity[df1$partisanship3 == "Republican"])
sd(df1$outlet_diversity[df1$partisanship3 == "Democrat"])
t_apa(
  res,
  es = "cohens_d",
  es_ci = FALSE,
  format = c("text"),
  info = FALSE,
  print = TRUE
)
df1 %>% cohens_d(outlet_diversity ~ partisanship3, var.equal = TRUE)

# slant diversity
res <- t.test(type_diversity ~ partisanship3, data = df1, var.equal = TRUE)
res
sd(df1$type_diversity[df1$partisanship3 == "Republican"])
sd(df1$type_diversity[df1$partisanship3 == "Democrat"])
t_apa(
  res,
  es = "cohens_d",
  es_ci = FALSE,
  format = c("text"),
  info = FALSE,
  print = TRUE
)
df1 %>% cohens_d(type_diversity ~ partisanship3, var.equal = TRUE)

# right extremity
res <- t.test(far_right ~ partisanship3, data = df1, var.equal = TRUE)
res
sd(df1$far_right[df1$partisanship3 == "Republican"], na.rm=TRUE)
sd(df1$far_right[df1$partisanship3 == "Democrat"],na.rm=TRUE)
t_apa(
  res,
  es = "cohens_d",
  es_ci = FALSE,
  format = c("text"),
  info = FALSE,
  print = TRUE
)
df1 %>% cohens_d(far_right ~ partisanship3, var.equal = TRUE)

# left extremity
res <- t.test(far_left ~ partisanship3, data = df1, var.equal = TRUE)
res
sd(df1$far_left[df1$partisanship3 == "Republican"], na.rm=TRUE)
sd(df1$far_left[df1$partisanship3 == "Democrat"],na.rm=TRUE)
t_apa(
  res,
  es = "cohens_d",
  es_ci = FALSE,
  format = c("text"),
  info = FALSE,
  print = TRUE
)
df1 %>% cohens_d(far_left ~ partisanship3, var.equal = TRUE)